How systems of single-molecule magnets magnetize at low temperatures 
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We model magnetization processes that take place through tunneling in crystals of single-molecule 
magnets, such as Mni2 and Feg. These processes take place when a field H is applied after quench- 
ing to very low temperatures. Magnetic dipolar interactions and spin flipping rules are essential 
ingredients of the model. The results obtained follow from Monte Carlo simulations and from the 
stochastic model we propose for dipole field diffusion. Correlations established before quenching are 
shown to later drive the magnetization process. We also show that in simple cubic lattices, m oc \/t 
at time t after H is applied, as observed in Fe8, but only for 1 + 21og 10 (/itj//i m ) time decades, where 
hd is some near-neighbor magnetic dipolar field and a spin reversal can occur only if the magnetic 
field acting on it is within some field window (— h w , h w ). However, the \fi behavior is not universal. 
For BCC and FCC lattices, m oat p , but p ~ 0.7. An expression for p in terms of lattice parameters 
is derived. At later times the magnetization levels off to a constant value. All these processes 
take place at approximately constant magnetic energy if the annealing energy e a is larger than the 
tunneling window's energy width (i.e., if e a > gfiBh w S). Thermal processes come in only later on 
to drive further magnetization growth. 

PACS numbers: 75.45.+j, 75.50.Xx 

Keywords: quantum tunneling, magnetization process, dipolar interactions, cooling history, dipole field 
diffusion 



I. INTRODUCTION 

Some magnetic clusters, such as Fes and Mni2, that 
make up the core of some organometallic molecules, be- 
have at low temperatures as large single spins. Ac- 
cordingly, these molecules have come to be known as 
single-molecule magnets (SMM's) A In crystals, magnetic 
anisotropy barriers, of energy U, inhibit magnetic relax- 
ation of SMM's, which can consequently proceed only by 
tunneling under the barriers. Magnetic quantum tun- 
neling (MQT) was first observed to take place through 
thermally excited states^ but temperature-independent 
"pure" MQT was observed shortly thereafter^ Dipo- 
lar interactions then play an essential role. They can 
give rise, upon tunneling, to Zeeman energy changes 
of nearly 1 K. This exceeds by many orders of mag- 
nitude the ground state tunnel splitting energies A 
that would follow for Fes and Mni2 from perturbations 
by anisotropics. Energy conservation would make pure 
MQT impossible for the vast majority of spins in the sys- 
tem. Hyperfine interactions between the tunneling elec- 
tronic spins of interest and nuclear spins open up a fairly 
large tunneling window of energy e w such that tunneling 
can occur if the Zeeman energy change 2sh upon tun- 
neling is not much larger than e w ~ More precisely, the 
tunneling rate T' for spins at very low temperature is 
given by 



r'(s h ) ~ r v(sh/s w ), 



(i) 



where T is some rate (whose value is not important for 
our purposes), rj(x) <~ 1 for | x |< 1, rj(x) ~ for x > 1, 
and e w > A. Other theories for MQT of SMM at very 
low temperatures have also been proposed^ We adopt 



Eq. Q here, regardless of theory or physical mechanism 
behind it. We let rj(x) = 1 for | x |< 1 and r)(x) = for 
x > 1, and refer to e w as the tunnel energy window. 

The interesting early time relaxation, m(0) — m{t) oc 
\fi, of the magnetization m of a system of SMM's 
that is fully polarized initially has been predicted^ ob- 
served experimentally^^ further explained^ and widely 
discussedA^ An unpredicted related phenomenon was 
later observed by Wernsdorfer et al.Ai the magnetiza- 
tion m of a system of Fes SMM's increases as y/i, where 
t is the time after a weak magnetic field is applied to an 
initially unpolarized system. There are important differ- 
ences between the early time relaxation of the magneti- 
zation of a system that is fully polarized initially^ and 
the magnetization process from an initially unpolarized 
stateiii Whereas the former effect depends crucially on 
system shape, the latter does not. Interesting questions 
arise: is this a universal effect to be found in all MQT ex- 
periments? If not, what does it depend on? How many 
time decades does the \/t regime cover? What is the 
final steady state magnetization? We have reported in 
Ref. [H a few results from MC simulations that answer 
some of these questions. The following problems were 
however not addressed in Ref. 12: (a) the crucial ef- 
fect that energy transfer (or the lack of it) between the 
magnetic system and the lattice has on the nature of the 
magnetization process; (b) a closely related problem, how 
the magnetization relaxes in zero field from an initially 
weakly polarized state. In addition, some results from 
our theory were given in Ref. IT2I but the theory itself 
was not. 

In this paper, we report extensive results from Monte 
Carlo (MC) simulations, an approximate theory, and re- 
sults that follow from it. The notation is first speci- 



fled. Unless otherwise stated, all energies and magnetic 
fields are given in terms of Ed and hd, respectively, where 
h c i = gfiBS/a 3 , a is the cubic lattice parameter, g is the 
gyromagnetic ratio, /is is the Bohr magneton, S is the 
spin size, and Ed — hdgp-BS. Temperatures are given in 
units of Ed/ks- For comparison purposes, the ordering 
temperature To and the ground state energy Eq are given 
in Table I for an Ising model with dipolar interactions 
in simple cubic (SC), body-centered-cubic (BCC), and 
face-centered-cubic (FCC) lattices. We let (h 2 } stand 
for the mean-squared value of the magnetic dipolar field 
h for random spin configurations, and a for l/[^/2np(0)], 
where p(0) is the probability density function (PDF) that 
the field on a randomly chosen site is when the spin 
configuration is random. 7 

The main results obtained follow. We show that the 
nonlinear in time magnetization arises from correlations 
that develop between spins and local magnetic dipolar 
fields, while cooling to low temperatures, before finally 
quenching to experiment. Furthermore, only the final 
energy —£ a reached just before quenching matters about 
the cooling protocol. More specifically, after quenching 
and applying a field H < 1 at t = 0, 

"»(*) " b TTK^ F (n, a/e w , a/ho), (2) 



where b 




for Tt < 1 

for 1 <Tt < (a/e w ) 1/p 
1 for (a/E w )^ p <Tt, 



(3) 
(4) 

(5) 



ho = (87r 2 /3 5 / 2 )gfj,BSpv, Pv is the number of spin sites 
per unit volume, p ~ 0.5 for SC lattices, and p ~ 0.7 for 
BCC and FCC lattices. These results, shown graphically 
in Fig. ^ are obtained from MC simulations in which 
£a > £w as well as from the theory given below. The 
theory also gives 



sin irp , — er 
i- = V2n—. 



P 



h 



(6) 



Additional results are mentioned in the plan of the paper 
below. 

The plan of the paper is as follows. The model and 
Monte Carlo simulations are described in Sect. [D] Wc 
also explain in Sect. [Q] how we simulate constant en- 
ergy processes. Section ITTT1 is devoted to the theory. As 
an introduction to the rest of the section, qualitative 
arguments are given in Sec. IIII Al The difference be- 
tween the probability density functions for up-spins and 
down-spins with a dipolar field h acting on them, which 
develop during the stage before quenching to very low 
temperatures, is derived in Sec. IIII Bl The theory for 
the magnetization process that takes place after a mag- 
netic field H is applied, soon after quenching, is given 
in Sec. IIII CI and the Appendix. Field functions that 
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FIG. 1: (a) ma(h 2 } /(e w e a H) versus Ft for an applied exter- 
nal field H — 1 and the shown values of e w . Symbols stand for 
Monte Carlo results, but dashed, full, and dot-dashed lines 
stand for results obtained from Eq. 121 H for e m — 0.05, 0.2, 
and 1, respectively. Before field H was applied, thermaliza- 
tion took place from an initially disordered configuration up 
to the time when the energy —e a reached the value —0.58. 
Data points for Ft < 4 follow from 4 x 10 4 runs for systems 
of 16 x 16 x 16 spins on a SC lattice. For Ft > 4, all data 
points are for averages over 4 x 10 s runs of systems of 8 x 8 x 8 
spins in a SC lattice, (b) Same as in (a) but for the shown 
values of e a and e w = 0.1. The continuous line is from the 
numerical solution of Eq. 12011 . (c) Same as in (a) and (b), 
but for different lattices. Symbols stand for MC results. Lines 
stand for results from theory: dashed is for SC and Fes lat- 
tices, dot-dashed and full lines are for BCC and FCC lattices 
respectively. For SC, BCC, FCC, and Fes lattices, partial 
thermalization previously took place at T a — 10, 20, 60, and 
1.06K, till e a = 0.31, 0.28, 0.53, and 0.15K, respectively. (For 
easy comparison with table I, units for Fes are given in Kelvin 
-but not for cubic lattices.) We used e w — 0.05 for all cubic 
lattices and e w = llmK for Fes. 



arc defined below develop in time holes that have been 
observed experimentallyii and are intimately related to 
the magnetization process. How these holes develop in 
time is the subject of Sec. IIII Dl In Sec. IIII El we study 
the relaxation in zero field of the magnetization of sys- 
tems that have been previously cooled in weak fields. In 
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TABLE I: vW)o stands for the rms spatial average of the 
dipolar field, a = [\/27rp(0)]~ 1 , andp(O) stands for the PDF at 
h — on a randomly chosen site, both for completely random 
spin configurations. Except for Fe8, they are both in units of 
gfj,BSa~ 3 , and follow from present MC simualtions; a is the 
lattice constant for cubic lattices. To and £o stand for the 
ordering temperature and ground state energy, respectively. 
Their values are taken from the shown references. ksT and 
£o are, except for Fes, in units of (gnsS) 2 a~ 3 . 



LATTICE 


n b 


V(h 2 )o 


a 


T Q 


£0 


SC a 


I 


3.655 


3.83(2) 


2.50(5) 


-2.68(1) 


BCC° 


2 


3.864 


4.03(2) 


5.8(2) 


-4.0(1) 


FCC a 


4 


8.303 


8.44(2) 


11.3(3) 


-7.5(1) 


Feg 6 


1 


46(1) mT 


31(1) mT 


0.4(1) K 


-0.51(2) K 



a These numbers for Tp and eq are from Re f. lUj 
' We assume an easy anisotropy axis as g iven in l3ll32l To 
and £o are from Ref. |30L 



Sec. IIII Fl we show how the magnetization crosses over 
to a linear in time behavior long after a field is applied if 
energy exchange between the magnetic system and heat 
reservoir takes place. Finally, the results obtained are 
discussed in Sec. IIVI 



II. MODEL AND SIMULATION 
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FIG. 2: Energy versus time in the annealing stage, when t < 
0, and after quenching, when t > 0. The time scale in the 
annealing stage, when t < 0, is such that the rate for energy 
lowering spin flips is 1. In the annealing stage T = 62 in both 
cases. For comparison purposes, the ordering temperature To 
and the ground state energy Eo are given in Table I. For t > 0, 
the continuous and dashed lines are for constant energy and 
for an isothermal process at T = 1, respectively. Equation 
(1) was enforced for both evolutions, but, as everywhere else, 
only when t > 0. The constant energy process was simulated 
assuming detailed balance with T u = 63. These data are for 
systems 16 x 16 x 16 spins, with e w = 0.1, averaged over 10 3 
runs. 



We use the MC method to simulate magnetic relax- 
ation of Ising systems of ±S spins, on simple cubic lat- 
tices with periodic boundary conditions (PBC), that in- 
teract through magnetic dipolar fields and flip under 
rules to be specified belowii^ In our PBC scheme, two 
spins interact whether they are in the same box or in 
different replicated boxes, but a spin at x,y,z interacts 
with a spin at x',y',z' only if —L x /2 < x — x 1 < L x /2, 
-Ly/2 < y-y' < L y /2, and -L z /2 < z - z' < L z /2, 
where L x , L y , and L z are the sides of the box-like sys- 
tems we simulate. Thus each spin interacts with N — 1 
spins, where N is the number of spins in the system. 
This scheme has been tested against a free-boundary 
scheme in which all spins in the system are allowed 
to interact, and found satisfactory when the system is 
only weakly polarized (or not at all) , as is expected from 
Griffith's theorem^ 

The system is first allowed to evolve towards thermal 
equilibrium at some "high" temperature T a . We assume 
ksT a ^ C//10, which implies that spin reversals then take 
place mostly through overbarrier processes. Accordingly, 
spin flips are then governed by detailed balance rules, and 
Eq. Q is not enforced. For reasons that become clear 
below, we also impose the restriction T a > To, where To is 
the long-range ordering temperature. One may think of 
this first process as a waiting stage that the systems may 
have to undergo in the cooling process, before quenching 
to a lower temperature where a tunneling experiment (as 
in Ref. Ill]) can be later performed. 



Let the time when this first stage is ended by quenching 
to a temperature below roughly 0.1L 7 / '(Sks) be t = 
Accordingly, Eq. is then enforced on all spin flips for 
t > 0. As for detailed balance, we then proceed as fol- 
lows. We may assume (1) that thermalization of a SMM 
system with the lattice does not take place (i.e., the en- 
ergy is constant) at very low temperatures, or (2) that 
heat is readily exchanged with the lattice, as is some- 
times the casein. In the latter case, we follow the stan- 
dard procedure. Monte Carlo results illustrate these two 
cases in Fig. [5] We fulfill the constant energy condition 
by enforcing detailed balance but using an appropriately 
chosen pseudotemperature T u . [From an expression be- 
low Eq. 0, k B T u w (h 2 )o/2e a . 2Q Note that T u > T a , 
since — e a cannot be smaller than the equilibrium energy 
at T a ]. We have checked that the mean energy is indeed 
constant under this rule, as illustrated in Fig. [21 Un- 
less otherwise stated (as in sect. IIII F|l . results reported 
below are for constant energy processes. 

MC results for the time evolution of 
m(h 2 )oa/ (e w e a H), after a field H = 1 is applied 
upon quenching, are shown in Fig. ^a) for various 
values of e w . Before quenching, the system was ther- 
malized at T a = lO/ks for some time till the energy per 
spin reached the value —0.58. Clearly, m scales with 
e w up to a crossover time of, roughly, T~ l a 2 / 'e^, where 
m levels off. Within the time range 1 < Ft < <7 2 /e 2 u , 
m oc \ft. Monte Carlo results that show how m scales 
with e a are exhibited in Fig. ^b) for e w = 0.1. Scaling 



holds also for the leveling off value m e of to. Equations 
(2-5) are inferred from these graphs, as well as from the 
theory given below, from which some of the constants 
are obtained. 



III. THEORY 

In this section, we try to understand the results of Sect. 
II and go somewhat beyond. 



A. Overview 

Rough arguments that can serve to introduce the 
derivations that follow in subsections below are given in 
this subsection. Before quenching to very low tempera- 
tures, the system is for some time at some temperature 
T a that is above the ordering temperature Tq. In addi- 
tion, ksT a > U /10. Overbarrier relaxation can then take 
place, and r = r exp(f7/fcsT) follows from Arrehnius' 
law. Consequently, spin flipping readily takes place in 
the laboratory within a second's time if To < 10~ 5 s. 
Some correlation between spin Sj and field hi at each site 
i can therefore be established, but no long-range order 
can develop as long as T a >T . Thus, immediately after 
quenching, the joint probability density to find a spin up 
Pl(h) with a field h acting on it is larger than Pi(h) if 
h > and vice versa. This is illustrated in Fig. Ola). As 
the system evolves towards equilibrium before quenching, 
we expect f(h) = P[{h) — p-\(h) to increase. It is shown 
in the following section and illustrated in Fig. E^b) with 
MC results, that thermalization effects prior to cooling 
to very low temperatures turn out to be rather simple: 
f(h) oc e a hp(h) if e a <| e \, where p(h) = p^{h) +Pi(h) 
and £o is the ground-state energy (—2.68 for a SC lat- 
tice). Since m = — J dhf(h) (from the definition of /), 
to oc e a is expected to ensue after quenching and appli- 
cation of an external field. 

After quenching in the laboratory to a temperature be- 
low approximately 0.1U / S, spin reversals take place only 
by tunneling between S z — —S and S z = +S states. For 
this, a spin must be within the "tunnel window" (TW), 
defined by rj in Eq. (1). Let some field H, fulfilling 
H 2 -c (h 2 )o, be applied after quenching. Then, h + H 
must be within the TW, that is, only spins with a mag- 
netic dipolar field within a TW centered on h = —H 
can flip. Consider some H > 0. A look at Fig. Ofa) 
shows that p±(—H) > p^(—H) then, which implies that 
to must increase with time after t = 0, since the number 
of upward flips minus the number of downward flips is 
proportional to p±(h) — P]{h) within the TW (the Zee- 
man energy for h = —H vanishes and we are assuming 
fcfiT ^> e w ). We therefore expect Pi(h) — ► p-\(h) with 
time if h is within the TW, but both pi(h) and p-\{h) to 
remain approximately constant for some time if h is out- 
side the TW. Monte Carlo results illustrate this behavior 
in Figs. |lfa) and (b). This is much as in the experi- 
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FIG. 3: (a) p a (h), versus h. All data points are for a system 
of 16 x 16 x 16 spins. □ and o stand for down-spin and 
up-spin, respectively. Full lines are for Gaussian functions 
centered on h = 0.84 and h — —0.84, with standard deviations 
of 3.83. The system was cooled to a temperature T a = 10, 
from an infinite temperature, and allowed to evolve towards 
equilibrium (but not to reach it) till its energy reached the 
value —e a = —0.42. For comparison purposes see To and 
the ground state energy values in Table I, and the energy is 
given by (1/2) J h\pf(h) — pi(h)]dh when no external field is 
applied. All points stand for averages over 1.5 x 10 5 runs, (b) 
\pi(h) — pi(h)]/(he a ) versus h for the shown values of e . No 
results are shown in the interval —0.1 < h < 0.1, as division 
by h in the h ~ neighborhood enhances errors greatly. 



ments of Wernsdorfer et aliii In addition, inspection of 
Fig. |3fa) suggests that f(h) oc h if h 2 <C (h 2 )o, whence 
to oc H also follows. Finally, since the number of spin 
flips increases linearly with e w if e 2 w <C (h 2 }o, we expect 
to oc £ w . This is all as in Eq. (2). 

The field h acting on a spin does not remain long within 
the TW after the spin flips. Changing fields produced by 
other spin flips bring h out of the TW, much as in a ran- 
dom walk. Thus, the "hole" (or "well") exhibited in Fig. 
IUa) widens as in a diffusion process. This is treated in 
detail in Sec. IIII Dl and illustrated in Fig. 0Jb). Initially, 
we expect the well's width 5w to increase linearly with 
time, since a dilute, spatially random distribution of spins 
gives a Lorentzian field distribution^ with a width that 
is proportional to the concentration. However, as the 
concentration c of flipped spins grows, the correspond- 
ing field distribution becomes Gaussian, with Sw oc ^fc. 
Consequently, we expect to oc \ft to ensue after some 
time since to oc — J dhf(h) and this integral grows as 6w 
after / nearly vanishes within the TW [see Fig. EJb)], 
and remains nearly constant therein subsequently. 

The to oc \ft growth stage comes to an end when Sw 
becomes as wide as the field distribution over all spins. 
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FIG. 4: (a) fa{h 2 )o/(£ a £wh), at the shown times, for a sys- 
tem of 16 x 16 x 16 spins, after it was thermalized at T = 10 
till e a — 0.51, then abruptly cooled, when H = 1 was ap- 
plied. Closed (open) symbols stand for MC (theory) results. 
A tunnel window s w =0.1 was enforced. All points stand for 
averages over 6 x 10 histories, (b) Same as in (a) but for 
the shown values of time. All results obtained from theory 
for f(h,t), that is, from Eq. 1201 were multiplied by 1.27 in 
parts (a) and (b). The need for this correction quite likely 
comes from second order effects in £ a /c not included in the 
theory (see Sec. I1II Bl . 

Then m levels off to m = m s . If however the magnetic 
system can exchange energy with the lattice at very low 
temperature, then a thermally driven magnetization pro- 
cess eventually takes over. This is the subject of Sec. 

mm 

B. Annealing 

Assume that either T a 3> To or that the time spent 
in the waiting stage is so short that the PDF that a 
randomly chosen site have field h is not drastically dif- 
ferent from the PDF, p[h), for a totally random spin 
configuration^ On the other hand, the conditional prob- 
ability to find ±5 on a site where the field is h ful- 
fills, in equilibrium, p(±S | h) oc exp(±ft,/fcsT). Now, 
since the joint distribution p(±S, h) for finding, on a 
randomly chosen site, h and ±S is in general given by 
P (±S,h)=p(±S\h)p(h), 

p{±S 1 h)otp{h)e ±h / kBT (7) 



follows in equilibrium. Monte Carlo results illustrate 
this point in Fig. EJa). Therefore, since the mean en- 
ergy is (E) = (1/2) J dhh[ Pl (h) - Pi(h)], it follows that 
(E) ~ (h 2 } /2k B T^. where (h 2 ) = Jdhh 2 p{h). The re- 
placement (/i 2 )o/2fcsT — * e a generalizes the above equa- 
tion to 

p(±S,h)ocp(h)e ±2hE « /{h2)o (8) 

for all times up to equilibration. Then, to leading order 
in e a h/(h 2 ) , 

hp 

Pl(h)- Pl (h)~2-^p(h), (9) 

where p-\{h) = p(S, h) and p±(h) = p{—S, h). Therefore, 
for Gaussian field distributions^ 

Mh )- Pl ( h )^JI^e' h2 ^\ (10) 
V 7T a 

since (/i 2 )o = c 2 then. All points for [p-f (h)— pi(h)]/(he a ) 
obtained from MC calculations for SC lattices collapse 
onto a single curve in Fig. OJb), in agreement with Eq. 

All results given below are for e a J a < 0.15. This can be 
accomplished by annealing at T > 4T . In addition, only 
applied fields much smaller than a are used. Thus, signif- 
icant higher-order (in e a h/ a 2 ) contributions are avoided. 

Note that Eqs. (|T|)- (|10f) are valid only in the annealing 
stage. They are inapplicable after quenching, because 
spins are then not free to adjust to local fields. 

C. Magnetization process 

We now examine the system's time evolution after 
cooling it abruptly, at time t — 0, to a temperature below, 
roughly, 0.1U/(Sk B )- A field H is applied for all i > 0. 
Then, real spin flips up to | S z |< S states can be ne- 
glected, and tunneling through the ground state doublet 
is the only available path for spin reversals. Accordingly, 
spin flips are allowed only if the spin's Zeeman energy is 
within the tunnel window. Now, if either the system is in 
thermal contact with a reservoir at temperature T such 
that k B T 3> e w , or the energy is constant and sufficiently 
high such that ksT u ^> s w , then 

m(t) = 2T J dhrj(H + h)f(h;t), (11) 

where f(h; t) = Pi(h; t) — p^{h; t). Let 

fi(h,H;t) = 2Tri(H + h)f(h;t), (12) 

whose physical meaning follows from comparison of these 
two equations. 

It is important to note that, Eq. I|12fl notwithstanding, 
fi =/= —f. This is because fi has to do with numbers of spin 
flips, not with changes in dipolar fields, which are brought 



about by such spin flips and also contribute to /. In order 
to establish how / evolves with time, it is convenient to 
define g(h;t) = f(h,Q) — f(h;t). Now, g(h;t) changes 
because spins which flip on sites where the field is within 
h and h+dh contribute ft(h, H; t)Atdh to gdh in time At, 
and also because of dipolar field changes brought about 
by such flips. These two effects are approximately taken 
into account in, 



dr / dhG{h,h;t-T)fi(h,H,T), (13) 



g{h-t) 



where G(h, h;t — r) is the PDF that, on a randomly cho- 
sen site, the field is h at time t, given that the field on 
that site was h at time r. Again, note that / = —ft 
would only follow if field redistributions, owing to spin 
flips, were disregarded, that is, if G(h,h;t — r) were re- 
placed by 5(h - h) in Eq. (JT3J. 

We next approximate G(h, h;t — r). Recall that H 2 <C 
(h 2 )o throughout. Consider a spin that points in opposite 
directions at times i and r. Let the fractional number of 
all such spins be 0(t, r). Since p(h,t) is approximately 
constant for the times of interest here, 0(t, r) ~ 0(t — t). 
For 0(f — t) <C I, we expecfe2i 



G(h,h;t-r) 



u(t 



?rp - h) 2 + u(t - r) 2 



(14) 



where u — 2ho<j), ho — (8ir 2 /3 5 > 2 )gpBSp Vl p v is the num- 
ber of spin sites per unit volume. Near the other end of 
the time range, when 0(t — r) — » 1/2, we lel^ 



G(h,h;t-r) 



-(h-h) 2 /(2v 2 ) 



2nv(t - t) 



(15) 



where v = cr^/20(t — r), when cj)(t — r) — > 1/2. 

To make progress, 0(i — r) must now be determined. 
Let p(h; t) be the PDF that a spin with a field h acting 
on it has flipped at least once in time t. In order to relate 
(f> and n, given by 



n(t) = dh p(h; t), 



(16) 



we note that the total number of spin flips is much larger 
than the number of first time flips at all times. This 
implies that spins flip several times before the dipolar 
field acting on them drifts away from the tunnel window. 
Accordingly, (see Fig. |SJ we adopt 



= n/2. 



(17) 



Finally, we write an equation for p(h; t), from which n(t) 
and, consequently, <f>(t) then follow. Arguing as for Eq. 



p(h;t)~J dr I dhG(h-h;t-T)v(h,H,T), (18) 
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FIG. 5: Fraction of all spins n that have flipped at least once, 
and quantity 0/10 versus Ft. In order to avoid cluttering, 
0/10, instead of (j> is plotted. All results are for e w = 0.1. 
o and □ stand for MC averages over 5 x 10 4 evolutions; the 
continuous and dashed lines follow numerically from Eq. (1201 . 
n/20 should equal 0/10, according to Eq. 117H . 



follows, where 

u(h, H, T ) = \p(h, 0) - p(h, T)] v (h + H), (19) 

and we have assumed once again that p(h; t) = p(h, 0). 
The desired equation for the magnetization process, 



dx 
~dl 



ci- 



p(-H) 
p(0) 



dx 



2e„ 



"' 2 / dT dru(H;t-r) 



•2e„ 



(20) 



where x — ma(h 2 )o/(£ a £wH) : c% = A\j2/ir, C2 = 2, and 
lu(H; t-r) = 1/G(-H, -H; t- r), is derived in the Ap- 
pendix. The above equation also holds for n(t) if we let 
x = na/e w , c\ = \/2/tt and c% = 1. 

For H 2 < (h 2 } and H < w(-H; t-r)^u{t~ 
t) in Eq. (JUJ and u(—H\ t - r) — > V2nv(t - r) in Eq. 
(|T5|l . Then, 



(21) 



is an approximation that fits data obtained from simu- 
lations of field distributions from various concentrations, 
20, of randomly placed spins, in an otherwise empty lat- 
tice, reasonably well. In Eq. (|2ip. u> depends on t — r 
through and, consequently, through n. 

The functional dependence of F, shown in Eq. (2), 
follows from Eqs and (JJTJ. To see how Eq. (3) 

comes about, notice first that u> <C e w when rt 1. 
Then, Eq. becomes dx/d(Tt) ~ c\ - c 2 x for H 2 < 
(h 2 ) a , whence Eq. (3) follows. For Eq. (4), we turn to 
numerical solutions of Eq. I|20|l which are exhibited in 
Figs. 1(a), (b), and (c). 

We can obtain an analytic expression for p from Eq. 
(|20p. Since numerical solutions of Eqs. (|20|l and l|21|) give 
x(t) oc t p for Ft > 1 and vanishingly small s w /a, we let 
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FIG. 6: Quantity p, defined by m oc t p , versus a /ho. Symbols 
stand for MC results for the shown lattices for e w — 0.05 and 
H < a/4. The lattice structure for Fes is as in Ref. |2^. The 
numbers next to symbols for the tetragonal lattice stand for 
the ratio of the basal to the z-axis lattice constant. (0,0,1) 
and (1,1,1) stand for easy magnetization axes. The full line 
stands for Eq. @. 



e.m/cr — > and try m oc t p as a solution for all t ^> T 1 . 
Equations JTJJ and (|2U give 



u> ■ 



irhgn 



(22) 



as e w /a — > 0. Consider first a; = na/e w , and, conse- 
quently, ci = \j2/n and C2 = 1. Then, for H <ti u, Eq. 
(|2*U1) becomes 



(7.T 



* efcr(r) 1 
ar — : . 

[i ar ax(t — t) + 1 ' 



(23) 



where a = who/ (2a). Assuming x ~ i p , the change of 
variable z = r/t, brings Eq. I|23l) . in the t — > oo limit, to 



o = W-- 

7T 



1 , pzP- 1 

dZ.—r- 
Oi{l 



(24) 



Equation |JBJ follows immediately from integral tables^. 
We now return to Eq. (|20() to work with x = 
ma(h 2 )o/(e a e w H), c\ — 4y%iVj anc ^ c 2 = 2. We now 
assume x ~ t p , and, proceeding as above, obtain (1) 
p = p and (2) m/n — > 2e a H/ (h 2 ) as i — > oo. Therefore, 
p in to oc i p is given by Eq. iffijl. 

Quantity p from Eq. ©, as well as results for p ob- 
tained from MC simulations for various lattices, are ex- 
hibited in Fig. El Monte Carlo generated data points 
are in reasonably good agreement with theory, that is, 
with Eq. ©, except for tetragonal lattices with 0.3 and 
0.4 ratios of the basal to the z-axis lattice constant (i.e., 
for a/ho — 1-08 and 1.42, respectively). This departure 
may follow from the fact that, for these rather asym- 
metric lattices, l/p(h — 0) turns out to be nonlinear 



in the concentration n of spin occupied lattice sites for 
n > 0.01. This is in contradiction with Lorentzian field 
distributions^! we have assumed for n <C 1 . Such distri- 
butions give p(h = 0) oc n _1 , which comes into Eq. 120|) 
through Eq. (JSTJ. 

The relation m/n — > 2e a H/ (h 2 )o holds just as well for 
small but finite e w /a, whence the constant (l/2)yJn/2 in 
Eq. (5) follows, since n — > 1 as t — > oo. In addition, an 
intuitive argument that gives some insight into Eq. (5) 
is given in Sect. IIIIFI 



D. Tunnel window's imprint and field diffusion 

It is interesting to see how g(h;t) behaves for early 
times, when u -C e w and <j> <C 1. It then follows from Eq. 
Ifnjl that G(h) can be replaced by S(h + H) in Eq. (fT3jl . 
which in turn yields 



g(h, H; t) = 2T V (h + H) [f(h, 0) - g(h; t)], 



(25) 



and, consequently, for Ft <C 1, 

g(h; t) = f(h, 0) v (h + H)(1- e- 2Tt ). (26) 

Therefore, when the variation of f(h, 0) over rj(h + H) is 
negligible, g(h;t) is the tunnel window's imprint. This is 
illustrated in Fig. |H[b). A conjecture to this effect was 
made in Ref. ITU In order to show how long the condition 
u <C e w , necessary for the validity of Eq. (|2*fj)l . holds, 
note that the total number of spins that flip in time t is 
approximately p(h; 0)e w Tt. It follows immediately that 
u < e w if Tt < a /ho. Therefore, Eq. holds while 

Tt < 13 This is illustrated in Fig. gfb). Monte Carlo 
results for exponential TW shapes that further illustrate 
this behavior are reported in Ref. |2JJ. 

After Ft ~ 1, g(h;t) no longer grows as in Eq. 126|l . 
with constant width. When Tt 3> 1, the width of g(h;t) 
grows with time. A spin flip at any one site induces 
variations of dipolar fields at every other site. Wherever 
a spin flipped at time r the field h evolves as in a random 
walk, away from h = —H at time t, and G(h, —H; t — r) 
is its distribution at time t, as given in Eq. (|14f) while 
4>(t — t) <C 1. In this sense, one may then speak of 
dipole field diffusion. It is the purpose of this section to 
calculate such field diffusion. A Fourier transformation 
of Eq. I|A1|I gives 



9(k;t) 



drG{k;t - r)n{bk) 



m(r) 
~~2f^' 

AC w 



(27) 



where g(k;t) — J dhexp(—ikh)g(h;t), G(k;t) — 

Jdhexp[-ik(h - bh)]G(h,h;t), n(bk) 
J dhexp[-ibk(h + H)]n{h + H), and b = 1 - 20(t - r). 
A similar equation [letting g(k;t) — > p(k;t) and 
to(t) — > in Eq. (|27|l ] obtains for p(k;t). Equation 
(12711 can be solved numerically with the help of Eqs. 
11141 l|21|) . Results obtained are shown in Figs. EJa) 
and (b). Perhaps it is worth pointing out that p(h) 
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where too is the initial magnetization, and similarly for 
Pl(h; 0) letting to — > —mo- For H <C 1, 



FIG. 7: (a) cr(l — m/mo)/e w versus Ft for zero field relaxation 
from weakly polarized states for the shown values of e w and 
of mo. Data points for Ft < 64 (for 64 < Ft < 4 x 10 4 ) 
follow from averages over 10 5 (some 2 x 10 4 ) MC simulations 
of systems of 16 x 16 x 16 (8 x 8 x 8) spins, (b) Same as in 
(a) but for BCC and FCC lattices, e w =0.1 and mo = 0.2. 
The points indicated by the horizontal arrows correspond to 
m — > 0. The sloping lines for t p , for the shown values of p, 
are given for comparison purposes. 



changes appreciably only after very long times when 
strong spin-spin correlations develop as the equilibrium 
ordered state is approached iS^L 

The spike in f/h, which develops at h = after some 
time, was not observed experimentally. 11 However, it is 
not a spurious effect. It is a consequence of the fact that 
the hole in / spreads out in time, which implies that 
/(0, t) < when t > 0. The spike in f/h at h = has a 
physical consequence: if H is switched off at some time 
t > 1, a hole in f{h,t) will develop at ft. = 0, and the 
magnetization will decrease in time accordingly^ 



E. Relaxation in zero field 

A slight variation of the experiment above is consid- 
ered in this section. Suppose a magnetic field H is ap- 
plied during the annealing stage, so that a polarization 
develops before quenching to very low temperatures, and 
let this applied field be switched off at t = 0, soon af- 
ter quenching. Note that this differs from the problem 
considered in Refs. and|t| of systems that are initially 
fully polarized. Now, for a weakly polarized system, we 
expect 



rn 



27R7 



(29) 



pt(M)-^(i 



m )p(h; 0), 



(28) 



follows. Now, formally, replacement of e a h/(h 2 ) by 
—mo/2 in Eq. © gives the equation above. Since every- 
thing else works for the relaxation of the magnetization 
the same as in Sec. IIII CI substitution of e a h/(h 2 )o by 
—mo/2 in Eq. (2) leads to 

m ~ to [1 — b'e w a F(Tt, cr/e w , a/ho)], (30) 

where b' = 2y/2/n and F is given by Eqs. (3)-(5). This is 
the desired equation. Consequently, to becomes vanish- 
ingly small after Tt ~ {a/e w ) 1 / p . As in sections above, 
p ~ 0.5 for SC lattices, but is different for other lattices. 
Up to then, m relaxes, after Ft i=s 1, as t p . Monte Carlo 
results illustrating this behavior are shown in Figs. EJa) 
and (b). 

Thus, systems that are weakly polarized at quenching 
relax afterwards in complete analogy to the magnetiza- 
tion process described in previous sections. Relaxation 
proceeds as t p , from Ft ~ 1 up to the time {a/e lu ) 1 / p , 
when to reaches a vanishingly small value. In addition, 
quantity p depends on lattice structure the same as above 
and as in Fig. |SJ This differs fundamentally from the 
relaxation of systems that are initially fully polarized, 
which are the subject of Refs. 5 and 7 8 9. In the latter 
condition: (1) the relaxation of the magnetization de- 
pends crucially on system shaped (2) a \fi relaxation is 
predicted to be independent of lattice structure; (3) this 
behavior holds while 5 ' 7 e w /a < Ft and to > 0.9to s , which 
take place much earlier than the magnetization processes 
that are the subject of this paper. 



F. Late time magnetization 

If the magnetic system does not exchange energy with 
the lattice, the magnetization finally levels off to the sta- 
tionary value m s ~ 1.6e a H/ (h 2 } , given by Eqs. (2) 
and (5). To gain some insight into this, note that: (1) 
to = J dhg(h); (2) g(h, t) comes near its stationary value, 
f(—H,0), in the vicinity of h = —H when Ft > 1; af- 
terwards, (3) J dhg(h,t) increases with time mostly be- 
cause the width 5w of g(h,t) increases up to the time 
when Sw ~ 2^J (h 2 )o, that is, when g(h,t) becomes as 
wide as p(h), and reaches its stationary state. Then, 
/ dhg(h) « f(-H, 0)2 V^o", which, making use of Eq. 
© for f(-H,0) with | H |< 1, and of a 2 - (h 2 ) , gives 
the desired estimate for m s . This comes to be when all 
spins have flipped at least once, i.e., when n(t) — > 1. The 
time evolution of n(t) is shown in Fig. |S]for e w = 0.1. 

On the other hand, even when UbT < 0.1U / S, heat ex- 
change between magnetic and lattice systems does take 
place in some systems of SMM'sJ^ Then, the magnetiza- 
tion can increase further with time as the system evolves 
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FIG. 8: m versus Tt for e w = 0.1, and H = 1. All data points 
are for MC simulations of 16 x 16 x 16 spin systems. Open 
(full) symbols stand for evolutions in which heat exchange 
between magnetic and lattice systems does not (does) take 
place. • is for e a = 0, A is for e a = 0.10, A is for e a = 0.097, 
□ is for e a = 0.29, ■ is for e a = 0.29, is for e a = 0.58, 
and ♦ is for e a = 0.58. Data points stand for averages over 
at least 10 3 MC runs. For comparison purposes, two straight 
dashed lines for m oc t v are shown. 

towards its thermal equilibrium state. As has previously 
been shown^ a relaxation time Tth can then be denned, 
which is given by 

r **~ (£)"■ (31) 

Accordingly, we expect 

H t 

m(t) ~ c (32) 

C T t h 

to hold when thermalization effects become dominant. 
Equation (|32[) . with c ~ 1.5, gives a rough fit of the 
nearly linear in time pieces in Fig. (jSJ). A crossover to 
this late time regime, which is driven by energy exchange 
with the lattice is easily appreciated in Fig. JHJ. The cor- 
responding crossover time t co is best expressed in terms 
of another crossover time, r s , when to oc y/i crosses over 
to m = m s . t co follows from Eqs. (l)-(5) and Eq. (|32|l . 
For e a ^> e w , the thermally driven linear behavior ap- 
pears after m levels off at m s . Letting a 2 ~ (h 2 )o and 
p»2, 



obtains. On the other hand, if e a <C s w , linear behav- 
ior appears early, and magnetization leveling off is pre- 
empted. By the same approximations, 




then follows. 



Even later yet, when t > Tth, the magnetization crosses 
over to another regime as long-range order sets in. This 
can be appreciated in Fig. 0for Tt > 10 4 . 



IV. CONCLUSIONS 

We have given MC and theoretical evidence to show 
that the to oc y/i behavior observed in experiments on 
Feg crystalsii after quenching and applying a small field 
H at t — is driven by correlations, which are previously 
established in the system while cooling to very low tem- 
peratures. Of the cooling protocol, only the correlation 
energy — e a at the time of quenching matters. The whole 
magnetization process is ruled by Eqs. (2)-(5) when the 
magnetic system is thermally isolated. 

The to oc y/i behavior has been shown to be nonuniver- 
sal. Both MC and theory show that this behavior ensues 
in SC lattices. Our MC simulations of an Fes model, 
in the appropriate lattice^ give m oc ^/i, in agreement 
with experiment. 11 However, MC simulations of FCC and 
BCC lattice systems give m oc t p , where p varies with lat- 
tice structure as illustrated in Fig. dc). Numerical solu- 
tions of Eq. (|20|l show that p depends on lattice structure 
through a/ho (ho, which is defined below Eq. (5), is pro- 
portional to the number of spin sites per unit volume) as 
shown in Fig. Values of a and p(0) are given in Table 
I for a few lattices. The m oc t p regime has been shown 
to cover the time range 1 < Tt < (a/ew) 1 ^. 

Similar results, namely Eq. (|30|l . have been obtained 
for the magnetization relaxation in zero field of systems 
that have been previously cooled in weak fields. This dif- 
fers from the relaxation of fully polarized initial states, 
where the 1 — to oc y/i behavior has been shown to be 
universal for 1 — to < 0.1. 5 More remarkably, our results 
differ qualitatively from the exponential relaxation pre- 
dicted in Ref. for weak field cooled systems. 

In order to test for system shape, size, and boundary 
effects, we have also performed MC simulations on spher- 
ical and box like systems of various sizes, with periodic 
and with free boundary conditions on SC lattices. Quan- 
tity p for box like systems with PBC appears to be size 
and shape independent as long as there are more than 
four spins on each side. Box like systems of L x x L x x L z 
spins with free boundary conditions give values of p that 
appear to agree in the L x — ► oo limit if L z > 16 as well as 
in the L z — > oo limit if L x > 8. Results for spheres with 
free boundary conditions also go to the same macroscopic 
limit : p = 0.49 ± 0.01 for SC lattices. 

The magnetization process depends significantly on 
whether the magnetic system is thermally isolated or not 
only after some crossover time r co , given by Eqs. t\'M and 
l)34p. This is illustrated in Fig. [H] It might be interesting 
to observe this effect on systems in which specific heat ex- 
periments have shown that energy is exchanged with the 
lattice within reasonable times at low temperatures^ 

The simple result of Sec. IIII Dl deserves a comment. It 
shows that experimentally observed holes in f(h,t) can 
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be used to establish the value of T, since holes develop 
independently of e w and of e a while Tt < 0.5, i.e., while 
/ at the bottom of the hole is approximately larger than 
one half its value at the top. This is useful because there 
is no other simple relation we know of that one can use 
in order to extract the value of T from othe experimental 
observations. For instance, application of Eq. (|2(jll to 
Fig. 4 in Ref. [ll] yields T" 1 ~ 24 s for Fe 8 . 

Except for the m oc \ft behavior, which has been ob- 
served in Fes over a limited time range^i, all the pre- 
dictions we make here have yet to be observed experi- 
mentally. Since time scales can be controlled, varying T, 
through the application of a transverse field, experimen- 
tal observation seems feasible. 
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APPENDIX A 



follows immediately, and, similarly for n(t) (letting m — > 
n, f — > p, and g — > p above). Substitution of Eq. I|A1(I . 
into Eq. I|A2|I . gives 



m ~ 2r 



2s w f(-H,0)- / dr-^W(t-r) 

JO 



(A3) 



where e w <C tr has been assumed, and 



W= J dh J dhr)(h + H)G(h,h;t - T)?](h + H). (A4) 



Approximations on W follow. First, T(t — t) <C 1 => 
u(t — t) <C e w , since u = 2h (f), n/2 < </> < n(l — n/2), 
and n(t) ~ (e w /a)Tt if Ft <C 1, whence G(h, h;t — r) — > 
5(h — h), and therefore 

W -> 2e w (A5) 
then. By the same argument, Y{t— r) 3> 1 u(t — r) 3> 
£„,, whence it follows that variations of G over the tunnel 
window are then negligible (from the assumption that 
£ w < a), and 



Equation Ijl^jl can be further simplified making use 
of the approximation fi(h,H,T) ~ + H)rh(r)/2s w . 
Then, 



m(r) 



dhG(h,h;t-r)r)(h + H) (Al) 



A similar equation follows for p(<), making the replace- 
ment v(h,H,T) -> r)(h + H)ti(t) /2e w in Eq. (I5)l. Now, 
from Eq. I|ll|) and the definition of g, 



Tfi — 2r / d/i[/(/i, 0) - g{h; t)]r)(h + H), 



(A2) 



W ^ (2e w ) 2 G{-H,-H;t-T) 



therefore. We now use 

W(t - t) 



2e„ 



[2e w G{-H,-H-t- T )]- 1 +1 



(A6) 



(A7) 



to interpolate between Eqs (|A5|) and JA6J. Substitut- 
ing this equation into Eq. (| A3|) [and into the analogous 
equation for n(t)] gives Eq. pO)l. which is the desired 
equation. 
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